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Abstract 

A new, "quantum mechanical" algorithm is proposed for automatic photopeak lo- 
cation in gamma-ray spectra from semiconductor and scintillator detectors. 



■o ', 1 Introduction 
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A great variety of programs can be found in the literature for automatic 
gamma-ray spectrum analysis with computers []T[. An important ingredient 
of every such program is some peak finding procedure. Usually the following 

<D two characteristics of a local maximum of continuous function are used in 

Ph! peak searching algorithms: 

^ ■ 1) Near the maximum the curve is convex and so its second derivative becomes 

negative, having a minimum value around the maximum of the peak. 
2) While passing a peak the first derivative changes sign. 

The use of the second derivative for automatic photopeak location was 
proposed by Mariscotti j2j. The real spectra, however, can have significant 
fluctuations because of the discrete and statistical nature of the data. So a 
discrete analog of the second derivative (the second difference) ddi = N i+ i — 
2Ni + Ni-i should be used, Ni being the pulse heights of the spectrum in 
channel i. 

To reduce further the influence of the statistical fluctuations, this second 
difference can be averaged and after n iterations replaced by the generalized 
second difference 

i+m jn+m h+m 

ddi(n,m)= E E ••• E 

j n =i-m j n -i=j n -m ji=32-m 



Its standard statistical deviation is sdi(n,m) = CjNj and therefore not 
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only ddi(n, m) < condition is demanded for the peak area location but also 
a big enough value of A H n ' m | 

° ° sdi{n,m) 

In |3j] some modification of this method was used. Instead of smoothed 
second difference one can consider the sum ddi(k) = Ejttfc CjNj where the 
Cj coefficients are chosen in such a way that ddi(k) < near the peak. For 
example in they take 
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p being some parameter to be optimized. Note that near the peak 
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In ||] the smoothed first derivative was used for photopeak searches. The 
peak is identified at the point where the first derivative changes sign from 
positive to negative. 

One can smooth not the first or second derivatives but the original data 
and establish the sign of the first derivative by simply comparing consecutive 
channels. To increase statistical reliability, one demands not only Ni > A^_i, 
but also N i+ i > Ni |5j] . As is claimed in , for low statistics the best result 
is given by the following smoothing 

Ni = l (N l+2 + 2JV<_i + 3N t + 2N. l+1 + N l+2 ) . 

One more method for photopeak position finding, although somewhat 
complex, is fitting the whole spectrum or its part by appropriate analytic 
function. A lot of various functions were suggested for this purpose, some of 
them modeling even Compton continuum jTj . If resolution is very good, as it 
is usually for semiconductor detectors, so that one has almost ideal gaussian 
peaks and the background in the peak area is nearly constant, it is possible to 
fit successive small sections of the data by Gaussian function and identify the 
peak position where the amplitude of the fitted Gaussian becomes significant 



In this note a new and simple enough algorithm is proposed for photopeak 
searches, which is completely different from the above described methods. 
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2 The idea 

The proposed method is based on a very simple but nice idea, which can 
be explained as follows. Suppose we place a small ball on the edge of the 
irregular potential wall: 



Classical ball will stop in front of the first obstacle. But if it is a quantum one 
and can penetrate through narrow barriers it still goes down to the potential 
wall bottom and oscillates there. 

For automatic photopeak finding this idea can be realized, for example, in 
such a way. Let us take a point somewhere on the right slope of the peak and 
let this point can jump left or right by one channel, or remain on its position, 



so that the probability to jump left is proportional to exp 
jump right exp 



-Ni 



-N 



and to 



N, 



j+r 



-Ni 



r- . This point will climb up quickly enough to the 

_WNi+i+Nij 

photopeak and will oscillate on the top. If now we follow up some amount 
of jumps and calculate the mean position, this just gives us the photopeak 
location. 

To smear the "quantum ball" even more, one can take for the probability 
to jump from channel i to channel i-1 : 



Pi,i-\ ~ E ex P 
k=i 



Ai 



N 



i—k 



where are some numerical constants. Changing them and the number m, 
one can govern the penetrating ability of the ball. 

3 Realization via discrete Markov chain 

In fact the above described set of channels and transition probabilities can 
be treated finite Markov chain @, [KJj . 
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For simplicity let us assume that from any state (channel) in this Markov 
chain only the closest left and right neighbour states can be reached with 
nonzero probabilities in one step, and that the probability for any state to 
remain unchanged is zero. Thus the transition probability matrix for this 
chain looks like 



P = 



( 

P21 






1 



-P32 





P23 








^34 









Where we shall take 



Pi,i±i = A J2 ex P 
k=i 



Ni 



±k 



1 



Ni 



(1) 



y/Ni±k + Ni 



(2) 



Ai normalization constant is defined from the Pi i-\ + P t 



1 condition. 



This Markov chain has a very simple invariant distribution [JTTJ : 

P12-P23 " ' Pn — 1,71 



Pl2 

U2 = ~p. — U\ 



PuP23 
U 3 = - - Ui 



Urn = 



^ 21 ^ 32 J 21 Mi,n-Kn-l,n-2 

Ui being defined from the normalization condition 

n 



p. 



-Ui 



21 



(3) 



(4) 



Now this invariant distribution has very sharp peaks which correspond 
to a local maximums in the original spectrum. This is illustrated by Fig.l. 
It is much simpler to identify peaks by some computer program in the u n 
distribution than in the original data. Fig. 2 gives a clear example of this 
(note the logarithmic scale on the vertical axis). 

In fact this method of photopeak finding can work even for very low statis- 
tics, as Fig. 3 shows. 



4 Photopeak fitting 

After the photopeak is located by the above described algorithm, one can fit 
it by gaussian function and extract its other characteristics. Fig. 2 and Fig. 3 
show that our " quantum" algorithm not only gives the peak position, but can 
also provide useful limits in which it is worthwhile to fit it by Gaussian. 
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For good enough gaussian peaks one can use a non-iterative method for 
the fast fitting (jTTj. Suppose we want to fit a Gaussian to the spectrum 
region from channel m\ to channel rri2- Let us form an array S(i) = In -J^ 

In J^' 1 , where m is some integer number. If we have a pure Gaussian 

1 *i+m — 1 



Ni = A exp 



2a 2 



then S(i) = doesn't depend on i. In fact this will not be the case because 
of statistical fluctuations and we should find a constant C which minimizes 
the sum 

1 1 



™g™ [S(i) - C} 2 

i=mi + l D 2 (i; 



1 1 



+ 



+ 



'7+777 Ni—i Ni +m -i 

where D(i) is a standard deviation for the casual quantity S(i). The solution 
of this linear minimization problem is well known 



C 



S(i) 



1712 — 171 

i=TOi+l 



m^—m 
^ 11- 



1 



(5) 



If we disregard statistical fluctuations of D(i), the statistical error in deter- 
mination of C is 

1712 — 171 1 1 _ ^l 2 

i=mi+l ^ \ l j 



So we get 



171 A 0" A ^ 



(6) 



The parameter m can be determined from the condition Act to be minimal. 
Let us take as an estimation D(i) = -l=. Then Act ~ — -, — = , which is 

v ' yNi m^/m%— mi— to' 

minimal for 

2 

m = -(m2 — mi) . (7) 



Analogously one can find the position and amplitude of the fitted peak by 

i-2p- 
2a 2 



considering arrays P(i) = ln^r = m(2t ^ +m) and Q(i) = IniV, = In A 



{l ^ . The results are 
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and 



lnA = 



TO2 In N- + ^ p) ' 



i=m\ 



1 



DN{i) 2 
'Aa\ 2 



V 



m 2 \ 
x 4 



a 



a 



(9) 



where m is defined by eq.7 and — <C relation was assumed while deriving 
the last equation. 

The degree of fit between the calculated and observed spectra can be 
further analyzed by means of some quantitative criteria, suggested in the 
literature []I2| . 



5 Conclusions 

It seems to us that the proposed algorithm is simple and effective enough to 
be recommended for practical applications. 

An earlier version of it was used with excellent results for uniformity stud- 
ies of about 1650 NaJ crystals for SND detector and also for vacuum 
phototriodes testing []T§ . 
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Figure captures 

Fig.l — Invariant probability distribution u(i) for single gaussian peak, m=3. 
Fig. 2 - - Invariant probability distribution for peak + background, m=3. 
Fig. 3 - - The same as on Fig. 2 for low statistics. 
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This figure "figl-l.png" is available in "png" format from: 



http://arXiv.org/ps/hep-ex/9506013vl 



This figure "figl-2.png" is available in "png" format from: 



http://arXiv.org/ps/hep-ex/9506013vl 



This figure "figl-3.png" is available in "png" format from: 



http://arXiv.org/ps/hep-ex/9506013vl 



